Dynamics of a Small Neutrally Buoyant Sphere in a Fluid and Targeting in Hamiltonian Systems 
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We show that, even in the most favorable case, the motion of a small spherical tracer suspended in a fluid 
of the same density may differ from the corresponding motion of an ideal passive particle. We demonstrate 
furthermore how its dynamics may be applied to target trajectories in Hamiltonian systems. 
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PACS numbers: 47.52.+j, 05.45.Gg, 45.20.Jj 

We show with the simplest model for the force acting on 
a small rigid neutrally buoyant spherical tracer particle in an 
incompressible two-dimensional fluid flow that tracer trajec- 
tories separate from fluid trajectories in those regions where 
the flow has hyperbolic stagnation points. A tracer will only 
evolve on fluid trajectories with Lyapunov exponents bounded 
by the value of its reciprocal Stokes number By making the 
Stokes number large enough, one can force a tracer in a flow 
with chaotic pathlines to settle on either the regular KAM-tori 
dominated regions or to selectively visit the chaotic regions 
with small Lyapunov exponents. These findings should be of 
interest for the interpretation of Lagrangian observations, for 
example in oceanography, and in laboratory fluid experiments 
that use small neutrally buoyant tracers. Moreover, since a 
two-dimensional incompressible flow is a particular instance 
of a generically chaotic Hamiltonian system, our results con- 
stitute a tool for targeting trajectories in Hamiltonian systems. 

Our starting point is the equation of motion for a small rigid 
spherical tracer in an incompressible fluid. 
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where v represents the velocity of the particle, u that of the 
fluid, Pp the density of the particle, pf, the density of the fluid 
it displaces, i^, the kinematic viscosity of the fluid, a, the ra- 
dius of the particle, and g, gravity. The terms on the right 
represent respectively the force exerted by the undisturbed 
flow on the particle, buoyancy. Stokes drag, the added mass, 
and the Basset-Boussinesq history force [[l|,||]. The terms in 
fl-^V^u are the Faxen corrections [|[|. The equation is as given 
by Maxey & Riley [Q], except for the added mass term, whose 
correct form was first derived by Taylor [^], as was pointed 
out by Auton et al. [^. The derivative Du/Dt is along the 
path of a fluid element 



Du du , 



(2) 



whereas the derivative du/dt is taken along the trajectory of 
the particle 



du du 

H ^ ~dt 



+ (v • V)u. 
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Equation (|I|) is valid where the particle radius and its Reynolds 
number are small, as are the velocity gradients around it, and 
the initial velocities of the particle and the fluid are equal. An 
excellent review of the history and physics of this problem is 
provided by Michaelides 

We shall simplify Eq. (|l[) even further with the aim of dis- 
covering whether in the most favorable case a tracer particle 
may always be faithful to a flow trajectory. With this in mind, 
we set Pp = Pf, so that the tracer be neutrally buoyant. At 
the same time we assume that it be sufficiently small that the 
Faxen corrections be negligible. Furthermore we exclude the 
Basset-Boussinesq term, since our approach — in which we 
follow Einstein and others [§-jl^] — is to obtain a minimal 
model with which we may perform a mathematical analysis 
of the problem. If in this model there appear differences be- 
tween particle and flow trajectories, with the inclusion of fur- 
ther terms these discrepancies will remain or may even be en- 
hanced |]rT[-pj|]. If we rescale space, time, and velocity by 
scale factors L,T = L/U, and U, we arrive at the expression 
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where St is the particle Stokes number St = 2a? U / {%vL) = 
2/9 (a/L)^Re/, Re/ being the fluid Reynolds number. The 
assumptions involved in deriving Eq. require that St <C 1. 

In the past it has been assumed that sufficiently small neu- 
trally buoyant particles have trivial dynamics [ p^ , p4] ], and the 
mathematical argument used to back this up is that if we make 
the approximation Du/Dt — du/dt, the problem becomes 
very simple: 
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u) = --St-i(v-u). 



(5) 



Thence v—u = (vo — Uq) exp(— 2/3 St ^ t), from which we 
infer that even if we release the particle with a different initial 
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velocity Vq to that of the fluid Uq, after a transient phase the 
particle velocity will match the fluid velocity, v = u, mean- 
ing that following this argument, a neutrally buoyant parti- 
cle should be an ideal tracer Although from the foregoing it 
would seem that neutrally buoyant particles represent a trivial 
limit to Eq. (jl|), this would be without taking into account that 
in a correct approach to the problem Du/Dt ^ du/dt. If we 
substitute the expressions for these derivatives given in Eqs. 
(|) and (|) into Eq. (Q), we obtain 

^ (v - u) = - ((v - u) . V) u - ^ St-i (v - u) . (6) 
We may then write A = v — u, whence 



dA 
~dt 



= -{J 



-St-i/) - A, 



(7) 



where J is the Jacobian matrix — we shall now concentrate 
on two-dimensional flows u — {u^, Uy) — 



J = 



If we diagonaUze the matrix we obtain 



dt 



A-2/3St-i 
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so if Re[X) > 2/3 St^ , Ad may grow exponentially. Now 
A satisfies det(J — AI) = 0, so A^ — tr J + dctJ = 0. Since 
the flow is incompressible, dxUx + dyUy = trJ = 0, thence 
— A^ = detJ. Given squared vorticity uo"^ = {dxUy — dyUx)^, 



and squared strain s 



si 



- where the normal component 



is si = dxUx~ dyUy and the shear component is S2 — dyUx 
dxUy, we may write 



= A2 = -dctJ= (s2-cj2)/4, 



(10) 



where Q is the Okubo-Weiss parameter [ p4|m ]. If Q > 0, 
A^ > 0, and A is real, so deformation dominates, as around 
hyperbolic points, whereas if Q < 0, A^ < 0, and A is imag- 
inary, so rotation dominates, as near elliptic points. Equation 
(0) together with dx/dt — A + u defines a dissipative dynam- 
ical system 



d^/dt = F(^) 



(11) 



with constant divergence V • F — 4/3 St^^ in the four di- 
mensional phase space ^ = {x,y, Ax, Ay), so that while small 
values of St allow for large values of the divergence, large val- 
ues of St force the divergence to be small. The Stokes number 
is the relaxation time of the particle back onto the fluid trajec- 
tories compared to the time scale of the flow — with larger St, 
the particle has more independence from the fluid flow. From 
Eq. (^, about areas of the flow near to hyperbolic stagnation 
points with Q > 4/9 St^^, particle and flow trajectories sep- 
arate exponentially. 
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FIG. 1. (a) Fluid trajectories (thick lines) and magnitude of Q 
(background shading: lighter is higher Q) for the time-independent 
model Eq. ( [l^ with A = 100 and B — 0. The flow is on a torus, 
(b) Separation of a neutrally buoyant particle trajectory (soUd line) 
with Stokes number St = 0.2 from the flow (dashed line) in regions 
of high Q allows the particle to wander between cells, (c) After a 
complicated excursion, a particle eventually settles in a zone of low 
Q. (d) The velocity difference — Ux between the particle and the 
flow against time for case (c). 

To illustrate the effects of St and Q on the dynamics of a 
neutrally buoyant particle, let us consider the simple incom- 
pressible two-dimensional model flow defined by the stream 
function 



'ilj{x,y,t) — Acos{x + B sinujt) cosy. 



(12) 



The equations of motion for an element of the fluid will then 
he X — dyip, y = —dxtp- tp has the role of a Hamiltonian 
for the dynamics of such an element, with x and y playing the 
parts of the conjugate coordinate and momentum pair. Let us 
first consider the simplest case where the time dependence is 
suppressed, by setting B ~ 0. Thence tp should be a con- 
stant of motion, which implies that real fluid elements follow 
trajectories that are level curves of ip. Such level curves are 
depicted in Fig. |l]a, which also shows contours of Q. The high 
values of Q are around the hyperbolic points, while negative 
Q coincides with the centers of vortices — elliptic points — in 
the flow. Figure |l]b shows the trajectory of a neutrally buoyant 
particle starting from a point on a fluid trajectory within the 
central vortex, with an inappreciable velocity mismatch with 
the flow. This mismatch is amplified in the vicinity of the hy- 
perbolic stagnation points where Q is larger than 4/9 St^^ to 
the extent that the particle leaves the central vortex for one of 
its neighbors, a trip that is not allowed to a fluid parcel. In the 
end a particle settles on a trajectory, proper for a fluid parcel, 
that does not visit regions of high Q. While this effect is al- 
ready seen in Fig. |l]b, it is more dramatically pictured in the 
trajectory shown in Fig. |l]c, in which the particle performs a 
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long and complicated excursion wandering between different 
vortices before it settles in a region of low Q of one of them. 
To illustrate the divergence of particle and fluid trajectories, 
and the fact that particle and fluid finally arrive at an accord, 
in Fig. |l]d we display the difference between the particle ve- 
locity and the fluid velocity at the site of the particle against 
time for the case of Fig. This difference is initially inap- 
preciable, and it converges to zero at long times, but during the 
interval in which the excursion takes place it fluctuates wildly. 
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FIG. 2. Poincare sections of trajectories in the time-dependent 
flowof Eq. withyl = 250, B = 0.3, and w = 1.0. (a) A chaotic 
fluid trajectory, (b) The motion of a neutrally buoyant particle with 
Stokes number St = 0.2 in the flow, (c) The velocity difference 
Vx — between the particle and the flow against time. 

Even more interesting is the case of time-dependent flow: 
£? 7^ in our model. As in a typical Hamiltonian system, as- 
sociated with the original hyperbolic stagnation points, there 
are regions of the phase space, which is here real space, dom- 
inated by chaotic trajectories. An individual trajectory of this 
kind, stroboscopically sampled at the frequency of the flow, is 
reproduced in Fig. Such a trajectory visits a large region of 
the space, which includes the original hyperbolic stagnation 
points and their vicinities where Q is large. Excluded from 
the reach of such a chaotic trajectory remain areas where the 
dynamics is regular; the so-called KAM tori. In our model 
these lie in the regions where Q < 4/9 St^'^. Now a neu- 
trally buoyant particle that tries to follow a chaotic flow path- 
line will eventually reach the highly hyperbolic regions of the 
flow. This makes likely its separation and departure from such 
a pathline, in search of another pathline to which to converge. 
However, convergence will only be achieved if the pathline 
never crosses areas of high Q. Figures^ and^ demonstrate 
this phenomenon: a particle was released in the chaotic zone 



with an inappreciable velocity mismatch. The particle fol- 
lowed the flow, until, coming upon a region of sufficiently 
high Q, it was thrown out of that flow pathline onto a long 
excursion that finally ended up in a regular region of the flow 
on a KAM torus. The regular regions of the flow then con- 
stitute attractors of the dissipative dynamical system Eq. ( pi] ) 
that describes the behavior of a neutrally buoyant particle. The 
chaotic trajectories in a Hamiltonian system are characterized 
by positive Lyapunov exponents. The Lyapunov exponents 
are an average along the trajectory of the local rate of conver- 
gence or divergence. Such a rate is measured by the quantity 
A. Hence, for a trajectory to be chaotic, it is a necessary con- 
dition that it visit regions of positive Q: an upper bound to Q 
is an upper bound to the Lyapunov exponent. 




FIG. 3. Small neutrally buoyant tracer particles with Stokes num- 
ber St — 0.2 collect in the centers of vortices in a two-dimensional 
turbulent flow simulation, (a) Q field at time t — 1 (lighter shad- 
ing is higher Q), (b) distribution at time t = 1 of particles that were 
uniformly distributed in the flow at time t — 0. 

We have considered the implications of these results for 
two-dimensional turbulent flows, in which Q defines three 
regions: in the vortex centers it is strongly negative; in the 
circulation cells that surround them, strongly positive, while 
in the background between vortices it fluctuates close to zero 
|[T^-|2C|]. We solve the two-dimensional vorticity equation 
for an incompressible fluid, dfOJ + J{ljj,'4>) — D^, where 
= —lyW^uj represents subgrid-scale dissipation; see, e.g., 
Mc Williams |21 1. We integrate this equation on a doubly pe- 
riodic domain using a pseudospectral method with 512^ col- 
location points and ly — 2.5 • 10^^; see Montabone [ p2| ] and 
Provenzale [ pO[ | for details. As a result of the dynamics, an 
initially uniform distribution of small neutrally buoyant par- 
ticles evolves in time towards an asymptotic distribution con- 
centrated in the inner part of vortices where Q < 0, and with 
voids in the areas crossed by fluid trajectories that visit re- 
gions where Q > 4/9 St^'^, as we illustrate in Fig. || This has 
importance consequences for the design of both fluid experi- 
ments using neutrally buoyant tracer particles, and Lagrangian 
drifters in geophysical flows. 

Now let us consider the following problem in the realm 
of Hamiltonian dynamics: given a generic chaotic Hamil- 
tonian H{pi, . . .pn,qi, ■ ■ ■ qn), we might like to find orbits 
with Lyapunov exponents with a given bound. In partic- 
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ular, we may be interested in locating small KAM tori in 
a chaotic sea. Inspired by the above results, we may fol- 
low the dynamics of a small neutrally buoyant hypersphere 
in a 2A^-dimensional hyperfluid. By this, we mean a parti- 
cle that follows a simplified Maxey-Riley equation, Eq. (Q), 
extended to 2N dimensions, which leads to a generalization 
of Eq. ([u]), dAo/dt = M ■ Ad, with elements Mij = 
Si.j [(— l)'^"'"A[(,j_|_i)/2] — 2/3 St^"'"], where [.] denotes the inte- 
ger part, replacing Eq. (^). Thus the Stokes number retains the 
same control over the Lyapunov exponents as in two dimen- 
sions. The quantity A may be viewed as the control signal: it 
vanishes when the desired trajectory is reached. In this light, 
the dynamics of a neutrally buoyant particle as demonstrated 
above may be thought of as a tool for targeting in Hamiltonian 
systems. 
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FIG. 4. (a) A particle in a hyperfluid defined by the 
Henon-Heiles Hamiltonian finds a KAM torus (arrowed) in the flow, 
(b) The velocity difference = Vx — Ux between the particle and 
the flow against time, (c) The Henon-Heiles energy E against time. 

In Fig. we illustrate our targeting mechanism with the 
Henon-Heiles Hamiltonian H = 1 /2{x^ + y"^ + + Py) + 
x^y — y^ /i. Figure ^ shows with a Poincare section how a 
particle released into a Henon-Heiles flow ends up on a KAM 
torus. In Fig. we plot one component of the velocity dif- 
ference between particle and flow; it displays two episodes of 
particle-flow separation before the particle settles on a KAM 



torus. Since it is not obliged to conserve the Henon-Heiles 
energy, the energy at which the particle finally settles into the 
Hamiltonian flow is a priori undefined; Fig. ^ shows a series 
of plateaux punctuated by rapid energy jumps correlated with 
the separations. If one is interested in finding a KAM torus 
at a given energy, this constraint should be imposed on the 
dynamics of the particle. 
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